Method of analyzing 3D geological structure using structure index

ABSTRACT

A method of analyzing a 3D geological structure using a structure index. The method includes the steps of estimating physical property values on common coordinates to calculate two or more physical property models on the same 3D grid (L×M×N); normalizing the physical property models, thus obtaining normalized physical property models which are then represented in a scatter plot of physical properties; converting distribution positions of the normalized physical property models on the scatter plot of physical properties into type angle (TA) and into type intensity (TI); determining a minimum of TI values, which classifies two or more classes on a scatter plot of the TA and TI values, to be a threshold; and analyzing the 3D geological structure based on a local extreme point and/or points in the scatter plot for the TA and TI using the threshold.

CROSS-REFERENCES TO RELATED APPLICATIONS

This patent application claims the benefit of priority from KoreanPatent Application No. 10-2011-0129545, filed on Dec. 6, 2011, thecontents of which are incorporated herein by reference in its entirety.

BACKGROUND OF THE INVENTION

1. Field of the Invention

The present invention relates to a method of analyzing a 3D geologicalstructure using a structure index (SI).

2. Description of the Related Art

Recently in order to acquire more pieces of information in survey areas,a variety of techniques are being carried out to conduct multi-surveys.This is because multi-surveys minimize the risk of uncertainty ofanalysis stemming from the use of a single survey technique and alsobecause the reliability of a final decision based on the survey resultsmay be increased. For this reason, thorough research into integratedanalysis to complement each data using multi-parametric geophysical dataand increase the reliability of analytic results is ongoing. However,conventional analysis into multi-surveys is merely confined to the levelof analyzing individual survey materials and then qualitativelycollecting analytic results, undesirably making it impossible toeffectively utilize the advantages of multi-surveys. This has beenconcluded to be because physical properties or resolutions measuredusing respective geophysical techniques are different from each other,and thus these are difficult to numerically summarize using a singledata processing technique. For example, many attempts have been made byresearchers (Vozoff and Jupp, 1975, Sill et al., 1977, Oldenburg, 1978)to make the sensitivity matrix (Jacobian matrix) of different surveydata into an inversion matrix so as to form joint inversion, but are notwidely utilized because solutions are numerically unstable and it isdifficult to make meaningful relationships.

In some embodiments of the present invention, however, integratedanalysis about the boundary information of geological structures, inlieu of a combination of physical properties, is carried out, thusovercoming the limitations of conventional integrated analysis. Such anintegrated analysis technique is advantageous because the acquired datais effectively used.

SUMMARY OF THE INVENTION

Accordingly, the present invention was developed keeping in mind theproblems occurring in the related art regarding multi-geophysical surveydata analysis, which led to the development of a method of effectivelyanalyzing multi-geophysical survey data so as to achieve undergroundgeological structure analysis corresponding to the basic purpose ofgeophysical surveys.

Typically when geological structures are analyzed using a singlephysical property, it is difficult to effectively carry out an analysisunder conditions of a small difference in the physical property betweenrock types, and thus some embodiments of the present invention provide amethod of more effectively analyzing a geological structure inconsideration of the correlation pattern of two or more physicalproperty models and the characteristics thereof.

In some embodiments of the present invention, an SI method which is ableto perform effective correlation analysis is provided.

The SI method consists of Type Angle (TA) and Type Intensity (TI)calculated using spatial correlations of two models, in which TAindicates the correlation pattern of two physical properties and TIindicates the distribution position in respective physical propertymodels.

The SI method for analyzing structures using TA and TI values may beeffectively employed when analyzing geological structures usingmulti-geophysical survey data. To further aid the understanding of theinvention, this method is illustratively applied to multi-geophysicalsurveys in the area of a volcanic caldera and the results thereof arerepresented.

An embodiment of the present invention provides a method of analyzing a3D geological structure using SI, including estimating physical propertyvalues (a density value and a resistivity value) on common coordinatesto calculate two or more physical property models (for example, adensity model and a resistivity model) D_(data)(L,M,N) andR_(data)(L,M,N) on the same 3D grid (L×M×N); normalizing the physicalproperty models, thus obtaining normalized physical property modelsND_(data)(l,m,n) and NR_(data)(l,m,n) which are then represented in ascatter plot of physical properties; converting distribution positionsof the ND_(data)(l,m,n) and NR_(data)(l,m,n) on the scatter plot ofphysical properties into TA using a four-quadrant inverse tangent andinto TI using a distance from an origin; determining a minimum of TIvalues, which classifies two or more classes on a scatter plot of the TAand TI values, to be a threshold; and analyzing the 3D geologicalstructure based on a local extreme point and/or points in the scatterplot for the TA and TI using the threshold. As such, normalizing may bedetermined by Equation 1 below:ND=[D−Mean(D)]/[(Max(D)−Min(D))/2]  <Equation 1>

wherein D is a physical property model D_(data)(L,M,N) orR_(data)(L,M,N), and ND is a normalized physical property model.

The method may further include comparison which is performed in such amanner that whether or not the two or more classes overlap with eachother is determined based on the determined threshold, betweendetermining the threshold and analyzing the 3D geological structure.

As such, when the two or more classes overlap with each other underconditions of being equal to or more than the threshold, determining thethreshold may be repeated so that the threshold is re-determined.

In an embodiment of the present invention, analyzing the 3D geologicalstructure may be performed by determining a rock type, and a structuredepending on physical properties and geological positions of theclasses.

The method may further include, before estimating the density value onthe common coordinates, surveying gravity at a position at which the 3Dgeological structure is to be analyzed, thus obtaining gravity data;correcting the gravity data, thus extracting a Bouguer anomaly; andperforming 3D density inversion using the Bouguer anomaly, thusobtaining the density model.

The method may further include, before estimating the resistivity valueon the common coordinates, surveying magnetotelluric (MT) at a positionat which the 3D geological structure is to be analyzed, thus obtainingMT data; correcting the MT data based on a remote reference, thusextracting MT data; and performing 3D MT inversion using the extractedMT data, thus obtaining the resistivity model.

In the present invention, the two or more physical property models onthe common coordinates may be estimated using results analyzed viainversion after the common coordinates have been set in the two or morephysical property models and by ordinary kriging.

Another embodiment of the present invention provides a method ofanalyzing a 3D geological structure using SI, including estimating afirst model value and a second model value on common coordinates tocalculate a first model D_(data)(L,M,N) and a second modelR_(data)(L,M,N) on the same 3D grid (L×M×N) having different physicalproperties, thus obtaining physical property models; normalizing thephysical property models, thus obtaining normalized physical propertymodels ND_(data)(l,m,n) and NR_(data)(l,m,n) which are then representedin a scatter plot of physical properties; converting distributionpositions of the ND_(data)(l,m,n) and NR_(data)(l,m,n) on the scatterplot of physical properties into TA using a four-quadrant inversetangent and into TI using a distance from an origin; determining aminimum of TI values, which classifies two or more classes on a scatterplot of the TA and TI values, to be a threshold; and analyzing the 3Dgeological structure based on a local extreme point and/or points in thescatter plot for the TA and TI using the threshold.

As such, normalizing may be determined by Equation 1 below:ND=[D−Mean(D)]/[(Max(D)−Min(D))/2]  <Equation 1>

wherein D is a physical property model D_(data)(L,M,N) orR_(data)(L,M,N), and ND is a normalized physical property model.

A further embodiment of the present invention provides a method ofanalyzing a 3D geological structure using SI, including estimating modelvalues on common coordinates to calculate two or more models on the same3D grid (L×M×N) having different physical properties, thus obtainingphysical property models; normalizing the physical property models, thusobtaining two or more normalized physical property models which are thenrepresented in a scatter plot of physical properties; convertingdistribution positions of the normalized physical property models on thescatter plot of physical properties into TA using a four-quadrantinverse tangent and into TI using a distance from an origin; determininga minimum of TI values, which classifies two or more classes on ascatter plot of the TA and TI values, to be a threshold; and analyzingthe 3D geological structure based on a local extreme point and/or pointsin the scatter plot for the TA and TI using the threshold.

As such, normalizing may be determined by Equation 1 below:ND=[D−Mean(D)]/[(Max(D)−Min(D))/2]  <Equation 1>

wherein D is a physical property model, and ND is a normalized physicalproperty model.

BRIEF DESCRIPTION OF THE DRAWINGS

The above and other objects, features and advantages of the presentinvention will be more clearly understood from the following detaileddescription taken in conjunction with the accompanying drawings, inwhich:

FIGS. 1A to 1C illustrate geologic map, 2-D gravity inversion result,and MT inversion result of the volcanic caldera used in an exampleaccording to an embodiment of the present invention;

FIG. 2 illustrates a terrain map and surveys points of the volcaniccaldera used in the example according to an embodiment of the presentinvention;

FIG. 3 illustrates Bouguer anomaly data used in the example anembodiment of the present invention, along with gravity and MTmeasurement points;

FIGS. 4A and 4B illustrate density models analyzed based on 3D densityinversion results according to an embodiment of the present invention;

FIGS. 5A and 5B illustrate estimates of the density model and theresistivity model in the same space according to an embodiment of thepresent invention;

FIGS. 6A and 6B illustrate explanatory diagrams of a process ofanalyzing a 3D geological structure using SI according to an embodimentof the present invention;

FIG. 7 illustrates a flowchart of the process of analyzing a 3Dgeological structure using SI according to an embodiment of the presentinvention;

FIGS. 8A and 8B illustrate normalized physical property models accordingto the present invention, in which FIG. 8A illustrates the scatter plotof physical properties and FIG. 8B illustrate TA and TI valuescalculated from the physical properties;

FIG. 9 illustrates 3D images of respective classes classified using theSI process according to an embodiment of the present invention;

FIGS. 10A and 10B illustrate the physical properties of respectiveclasses classified using the SI process according to an embodiment ofthe present invention;

FIGS. 11A and 11B illustrate the distributions of physical properties ofrespective classes classified using the SI process according to anembodiment of the present invention;

FIG. 12 illustrates a flowchart of a process of analyzing a 3Dgeological structure using SI according to an embodiment of the presentinvention; and

FIG. 13 illustrates a flowchart of a process of analyzing a 3Dgeological structure using SI according to an embodiment of the presentinvention.

DESCRIPTION OF SPECIFIC EMBODIMENTS

Hereinafter, embodiments of the present invention will be more fullyunderstood with reference to the appended drawings. Throughout thedrawings, the same reference numerals are used to refer to the same orsimilar elements. Furthermore, descriptions of known techniques, even ifthey are pertinent to some embodiments of the present invention, areregarded as unnecessary and may be omitted when they would make thecharacteristics of an embodiment of the invention and the descriptionunclear.

FIGS. 1A to 1C illustrate geologic map, 2-D gravity inversion result,and MT inversion result of the volcanic caldera according to anembodiment of the present invention, and FIG. 2 illustrates a terrainmap thereof.

As illustrated in FIGS. 1A to 1C and 2, the volcanic caldera used in anexample is the oval-shaped ring vent caldera positioned at the Euisungsub-basin in Korea, with a size of 16(EW)(N−S) km.

FIG. 1B illustrates 2D gravity survey results, and FIG. 1C illustrates2D MT survey results performed in the same area.

The volcanic caldera used in an example is known to have been formed viacentral eruption, sedimentation, formation of ring faults, and eruptionalong the faults in the tertiary cenozico era of about 54.5 Ma (YoonSung-Hyo, 1988), but the accurate production mechanism or lowerstructure thereof has not yet been ascertained. The sedimentary rockswithin the volcanic caldera are weakly folded and provided in the formof hornfels, and the central volcanic complexes which are rhyolitic areplaced on the Hayang group. Also, the sedimentation pattern is seen inthe central volcanic complexes, which means that additional internaleruptions occurred in the central portion and accompanied by hypabyssalmagma tism and sedimentation, resulting in fault formation (YoonSung-Hyo, 1988; Jang Gi-Hong and Park Sun-Ok, 1997). The overallstratigraphic structure of the volcanic caldera area includesNagdong-Hasandong-Jinju formations (Sindong Group), andIljik-Hupyeongdong-Jeomgog-Sagog-Banyaweol-Hwasan-Chunsan formations(Hayang Group), which are sequentially positioned upwards, and thecentral rhyolitic volcanic complexes and igneous intrusions are known tobe formed finally (Park Gyesoon et al., 2008).

Below is a description of a geophysical survey data conducted to analyzea 3D geological structure according to an embodiment of the presentinvention.

To analyze the 3D geological structure of the volcanic caldera locatedat the Euisung sub-basin in Korea, gravity survey (S10) and MT survey(S20) are first performed.

The measurement points for gravity survey are represented as black dotson the to terrain map of FIG. 2, and survey was performed at 510 pointsat intervals of about 1 km in the range of about 32 km around thevolcanic caldera.

The gravity data measured upon gravity survey (S10) is then subjected togravity correction (S11).

As such, gravity correction (S11) enables Bouguer anomaly data to beextracted.

To conduct gravity correction (S11), a typical gravity correctionprocedure, including latitude correction, tide correction, driftcorrection, free-air correction, Bouguer correction and terraincorrection may be performed, thereby extracting gravity effects due tounderground density structures.

FIG. 3 illustrates Bouguer anomaly data extracted via gravity correctionaccording to an embodiment of the present invention, together withgravity and MT measurement points.

As illustrated in FIG. 3, the Bouguer anomaly data is low in thenortheastern and southwestern parts of the survey area. The northeasternpart has a small number of gravity measurement points and thus theresulting estimation error is considered to be partially contained, andthe low anomaly of the southwestern part is analyzed to be due to thelow density of Palgongsan granite. Also, the volcanic caldera shows ahigh gravity anomaly zone at the center thereof, which is considered tobe due to the effects of igneous rocks intruded along the ring fault,and the low anomaly zone of the middle of the caldera is analyzed to bedue to the structure resulting from sedimentation of pyroclasticsedimentary rocks.

The Bouguer anomaly data extracted via gravity correction (S11) is thensubjected to 3D density inversion (S12) using a Marquardt-Levenbergmethod, thus extracting a density model.

In order to stably perform 3D density inversion (S12), a restrictionmatrix is added to the inversion so as to obtain similar densitydistributions between adjacent blocks, and the lower inversion block isset to be greater than the upper inversion block. Furthermore, the rangein which the density changes is set to be 0.2˜0.3 g/cm (Park Gyesoon etal., 2008).

FIGS. 4A and 4B illustrate the density models analyzed based on the 3Ddensity inversion results according to an embodiment of the presentinvention.

The density models analyzed based on the inversion results are depictedin FIGS. 4A and 4B. FIG. 4B illustrates the 3D image of the distributionof a low-density body, which shows the southwestern Palgongsan graniteand the pyroclastic sedimentary rocks of the middle of the caldera, andthe northeastern low-density body is analyzed to be a structure which ismerely numerically represented because there are insufficientmeasurement points.

MT survey (S20) is performed under conditions of the followingmeasurement points and system.

In the example according to an embodiment of the present invention, MTsurvey was performed at 32 points represented by the triangular shape inFIGS. 2 and 3 using a MTU-5A system. The distance between themeasurement points is about 2˜3 km, and the measurement points areradially distributed around the volcanic caldera. In the MT survey, theelectric fields in two directions of N-S and E-W and the magnetic fieldsin three N-S, E-W and vertical directions were measured in the bandrange of 0.001-320 Hz. To minimize noise effects, measurements weretaken for 15 hr during the nighttime when electromagnetic noise is low.

The MT data extracted via MT survey (S20) is then subjected to MTcorrection (S21) based on a remote reference.

Upon MT correction (S21), data obtained from the continuously operatingreference station (CORS) of Esashi, Japan which is about 1,000 km awayfrom the volcanic caldera area was used as the remote reference.

Furthermore, the MT data is corrected using an MT-Editor program whichis provided together with the MTU-5A system.

The MT data corrected via MT correction (S21) is then subjected to 3Dinversion (S22), thus extracting a resistivity model.

The 3D MT inversion according to some embodiments of the presentinvention adopts WSINV3DMT code (Siripunvaraporn et al., 2005) which isexpanded from a 2D inversion method, namely, Occam's inversion(Siripunvaraporn and Egbert 2000).

To analyze the geological structure using the spatial correlationbetween the density model extracted via 3D density inversion (S12) andthe resistivity model extracted via 3D MT inversion (S22), two modelshaving the same coordinates are required.

To this end, common coordinates are set and a density value and aresistivity value are respectively estimated on the common coordinatesin S100 and S200.

In the example of an embodiment of the present invention, estimation ofthe density and resistivity values on the common coordinates may beconducted using ordinary kriging based on S-GeMS (StanfordGeostatistical Modeling Software).

FIGS. 5A and 5B illustrate the density model and the resistivity modelestimated in the same space according to an embodiment of the presentinvention.

FIGS. 6A and 6B illustrate explanatory diagrams of the process ofanalyzing the 3D geological structure using SI according to the presentinvention.

As illustrated in FIGS. 6A and 6B, the method of analyzing the 3Dgeological structure using SI according to an embodiment of the presentinvention classifies structures using TA calculated from the correlationof two physical property values in the same space via the normalizedphysical property models and TI determined by the distribution positionof the corresponding value in the same physical property model.

A projection technique such as UTM (Universal Transverse Mercator) isused to classify structures, so that TA and TI are represented just aslatitude and longitude values and thus the structures are classifiedalong curves shown on the basis of the local extream points.

Basic suppositions proposed to use the SI method according to anembodiment of the present invention are described below.

1) The physical properties analyzed in respective physical propertymodels are different from each other but show the same undergroundgeological structure. As a result of geophysical surveys, geologicalstructures having different physical properties are separately analyzed.

2) The different physical property models change with the spatialcorrelation depending on changes in the geological structure. Becausethe different physical property models show the same undergroundgeological structure based on the above 1), they have a spatialcorrelation.

3) The different physical property models have similar survey depths andspatial resolutions. Because the structure is analyzed using thecorrelation of physical properties, the different physical propertymodels should share a similar survey space.

4) In respective physical property models, the anomaly value equal to ormore than a specific value shows a specific geological structure. Evenwhen the geological structures are different from each other, they mayhave similar physical property values. In this case, it is verydifficult to analyze the structure using geophysical surveys. Thus insome embodiments of the present invention, respective physical propertymodels are normalized, and if any model does not have the specificanomaly value as in the structure No. 0 of FIG. 6B, it is regarded as aboundary or variation zone of individual geological structures and thusis not classified into a specific geological structure and the analysisthereof is excluded.

Based on the above suppositions, the analysis of the geologicalstructure using SI is described below.

FIG. 7 is a flowchart illustrating the process of analyzing the 3Dgeological structure using SI according to an embodiment of the presentinvention.

Some embodiments of the present invention adopt an analysis method usingtwo physical properties, for example, density and resistivity.

As illustrated in FIG. 7, the density model extracted via 3D densityinversion (S12) and the resistivity model extracted via 3D MT inversion(S22) are subjected to density value estimation (S100) and resistivityvalue estimation (S200) on the common coordinates, respectively, usingkriging which estimates data values at any position in the context ofthe spatial distribution, so as to calculate the density modelD_(data)(L,M,N) and the resistivity model R_(data)(L,M,N) on the same 3Dgrid (L×M×N), thereby acquiring the desired data.

Subsequently, normalization (S300) is performed, and specifically,taking into consideration the distribution of respective physicalproperties, the density model D_(data) (L,M,N) and the resistivity modelR_(data) (L,M,N) present on the same grid are normalized using thefollowing Equation 1, and ND_(data) (l,m,n) and NR_(data) (l,m,n) areplotted in a scatter plot of physical properties.ND=[D−Mean(D)]/[(Max(D)−Min(D))/2]  <Equation 1>

wherein D is the physical property model D_(data) (L,M,N) or R_(data)(L,M,N), and ND is the normalized physical property model.

After normalization (S300), TA and TI conversion (S400) is conducted insuch a manner that respective positions of ND_(data) (l,m,n) andNR_(data) (l,m,n) distributed on the scatter plot of physical propertiesare converted into TA using four-quadrant inverse tangent and into TIusing the distance from the origin.

Subsequently, threshold determination (S500) is performed. Specifically,as illustrated in FIGS. 6B and 8B, the effects due to the same structurewhich may be estimated from the physical property values on the scatterplot of the TA and TI values enables the formation of an independentsingle class. Based thereon, the minimum of the TI values, which mayclassify two or more classes, is determined to be the threshold.

Subsequently, comparison (S600) is performed in such a manner that,based on the threshold determined in S500, whether the classes areclearly classified on the scatter plot of the TA and TI values underconditions of being equal to or more than the threshold is determined.

In the comparison (S600), the case where the classes are unclearlyclassified under conditions of being equal to or more than the thresholdturns back to the threshold determination (S500) so that the thresholdis re-determined.

When the classes are clearly classified under conditions of being equalto or more than the threshold in the comparison (S600), classdetermination (S700) is performed, specifically the corresponding valueis finally decided as the threshold and the classes are fixed on thescatter plot of the TA and TI values.

Among the classes acquired in S700, respective parabolas havingindependent maxima in the TI range equal to or above the threshold areclassified based on TA, and structure analysis thereof is carried out.

Specifically, 3D geological structure analysis (S800) is performed insuch a manner that the rock types and structures are determined inconsideration of the physical properties and the geological positions ofthe classes, and the spatial distribution of respective classes areformed in 3D to analyze the spatial characteristics of the geologicalstructures.

Example

In order to apply the method of analyzing the 3D geological structureusing SI according to the present invention, the density and resistivitymodels of the volcanic caldera distributed on the same coordinatesillustrated in FIGS. 5A and 5B are used.

The normalized physical property models are represented in a scatterplot of physical properties as depicted in FIG. 8A, from which TA and TIare calculated as illustrated in FIG. 8B.

To ascertain the applicability of the SI method according to anembodiment of the present invention, geological classes are classifiedat all the minimum values as possible as seen in FIG. 8B, and toeffectively analyze the geological structure of the caldera, thethreshold is set to 0.4 at which the parabolas are separated. The rangesof physical properties of the classes are shown in Table 1 below.

TABLE 1 Physical Properties of Respective Classes by SI NormalizedResistivity Normalized relative Structure Index Resistivity (log(Ω · m))Density Relative Density (g/cm³) Class 1 −0.5366~+0.5355 0.99~2.80−0.9882~−0.2023 −0.043~+0.076 Class 2 −0.6627~−0.2841 2.06~2.79−0.5236~−0.2034 −0.057~−0.015 Class 3 −0.7607~−0.3482 2.41~3.02−0.3718~−0.1064 −0.070~−0.022 Class 4 −1.1862~−0.3885 2.80~3.63−0.2000~+0.1624 −0.115~−0.027 Class 5 −0.6582~−0.2613 3.52~4.16+0.1130~+0.3911 −0.057~−0.013 Class 6 −0.3176~−0.0722 3.98~4.55+0.3135~+0.5623 −0.019~+0.008 Class 7 −0.1018~+0.8138 3.43~5.58+0.0729~+1.0118 +0.005~+0.011 Class 8 +0.3493~+0.7074 2.58~3.47−0.2962~+0.0931 +0.055~+0.095

3D images of respective class structures classified using the SI methodare depicted in FIG. 9. Considering prior research materials (YoonSung-Hyo, 1988, Park Gyesoon et al., 2008, Yang Jun-Mo et al., 2008) andspatial anomaly areas of physical property models, the volcanic calderaarea may be largely classified into three types of geological structuresvia a geophysical survey, as illustrated in FIGS. 9 and 10A and 10B.

1) Classes 2, 3, and 4: Palgongsan Granite and Pyroclastic SedimentaryRocks above the Center of the Volcanic Caldera

This geological structure has a low density which may decrease up toabout 0.1 g/cm³ and a low resistivity of about 200 ohm-m. This structureis shown near the center of the volcanic caldera which extends to adepth of about 1 km from the surface of the earth, which includes thepyroclastic sedimentary rocks and the Palgongsan granite. From the factthat the pyroclastic sedimentary rocks are observed up to a depth ofabout 1 km, it is concluded that the sedimentary rocks were formed bythe collapse of the middle of the volcanic body after an initialeruption. However, because the southwestern part of the survey area inwhich Palgongsan granite is observed not to have any MT survey points,it is difficult to have low resistivity as in the analytic results, andon the contrary, this structure is regarded as having high resistivitybecause of the properties of granite which has high resistivity. If theMT survey is carried out near the Palgongsan granitic area, thePalgongsan granite and the pyroclastic sedimentary rocks will be able tobe analyzed and shown to be different structures.

2) Class 6: Igneous Rocks Intruded in the Ring Fault of the VolcanicCaldera and Intrusive Granite over the Northeastern Part of the Caldera

This geological structure has high density and resistivity exceeding 0.1g/cm³ and 10,000 ohm-m. Such physical properties show the properties ofigneous rocks of this area, and manifest the intrusive igneous structureover the ring fault of the volcanic caldera and the northeastern part ofthe survey area.

3) Class 7: Basement Rocks

This geological structure has high density and comparatively lowresistivity. Considering the spatial position of this structure,physical property values, prior research materials, and the geologicalstructure around the volcanic caldera, the above structure is analyzedto be basement rocks or Sindong group of FIGS. 1B and 1C, and has a highdensity but comparatively low resistivity with respect to the geologicformation having superior electrical conductivity at the lower portionof the Euisung sub-basin proposed by Lee C. K (2006). If analysis ismerely based on respective physical property models, this class is noteasy to classify apart from class 6 having high density and highresistivity, but the use of the SI method facilitates the classificationthereof.

However, as shown in FIGS. 8B and 10A and 10B, class 5 is analyzed to bean intermediate zone of class 4 and class 6 and is not included instructure analysis. Also, the geological structure of class 1 may bedifferent from the others, but is shown at the edge of the survey areadue to the absence of gravity and MT measurement points, and thus datafor the area in which the reliability of inversion is low may beeffectively removed using the SI method.

FIG. 11A illustrates the density and resistivity distributions of theentire survey area and respective classes, and FIG. 11B illustrates theprobability distributions of the density and resistivity of the entiresurvey area and respective classes. From this, the range of physicalproperty values of respective classes and the number of outcropsdepending on the physical properties may be checked, and thus thephysical properties of respective classes may be ascertained at firstglance. Also, in the probability distributions of the entire surveyarea, the physical property values are dispersed over a wide range,whereas the structures of the classes are divided and show the physicalproperty values in narrow ranges. Furthermore upon structure analysis,there may be seen effects in which structures are similar in terms of asingle physical property value but may be different in terms of anintegrated analysis of two physical properties.

FIG. 12 illustrates a flowchart of a process of analyzing a 3Dgeological structure using SI according to another embodiment of thepresent invention.

In this embodiment of the present invention gravity data and MT data areused as different physical properties. Even when different physicalproperties other than the gravity data and the MT data are utilized,they may be included in the scope of the present invention.

FIG. 13 illustrates a flowchart of a process of analyzing a 3Dgeological structure using SI according to a further embodiment of thepresent invention.

This embodiment of the present invention illustrates gravity data and MTdata as different physical properties. Even when different physicalproperties other than the gravity data and the MT data are utilized,they may be included in the scope of the present invention.

As described hereinbefore, the present invention provides a method ofanalyzing a 3D geological structure using SI. According to the presentinvention, the SI method consists of TA and TI values calculated usingthe spatial correlation of two models in which TA indicates thecorrelation pattern of two physical properties and TI indicates thedistribution position in respective physical property models. The SImethod for analyzing the structure using these TA and TI values enablesthe effective analysis of geological structures using multi-geophysicalsurvey data.

When this method is illustratively applied to imaging the 3D geologicalstructure of the volcanic caldera, effective analysis can be found to bepossible.

Also because the SI method is used to classify structures using TA andTI values calculated from the spatial correlation of different physicalproperties and the distribution thereof, any structure, which isdifficult to classify due to similar physical properties when using anindependent geophysical survey analysis method, may be classified andanalyzed using different physical properties, and the SI method can beeffectively employed in 3D imaging.

The 3D underground structure of the volcanic caldera analyzed usingdifferent physical properties, for example gravity data and MT data,according to an exemplary embodiment of the present invention correlateswell with the wide underground structure resulting from surface geology,existing gravity and magnetic surveys. First, the spatial distributionof intrusive igneous rocks with high density and resistivity valuesdisposed at the ring fault of the volcanic caldera and in thenortheastern part of the survey area, and of the pyroclastic sedimentaryrocks above the middle of the caldera, which extends to a depth of about1 km, was analyzed in 3D. Also, the lower basement rocks having acomparatively low resistivity and high density were observed near about5 km. These results are in good correlation with those of previousresearch papers.

The SI method according to some embodiments of the present invention canobjectively analyze changes in physical properties of each structure andthe distribution thereof, and can extract useful information which wasnot found via individual analysis, and these results are able to serveas basic data for the precise re-analysis of local areas.

Moreover, because the SI method is determined using vectors, it is easyto extend it.

Therefore, even when two or more multi-geophysical surveys are carriedout, an integrated analysis can be effectively achieved, and a varietyof geophysical survey data in the same area in the future can beobjectively and effectively stored in the form of a database and can beeffectively applied to the imaging of underground models.

Although the preferred embodiments of the present invention have beendisclosed for illustrative purposes, those skilled in the art willappreciate that various modifications, additions and substitutions arepossible, without departing from the scope and spirit of the inventionas disclosed in the accompanying claims.

What is claimed is:
 1. A method of analyzing a three-dimensional (3D)geological structure, comprising: estimating a density value and aresistivity value on common coordinates L, M, and N to calculate adensity model D_(data) (L,M,N) and a resistivity model R_(data) (L,M,N)on a same 3D grid (L×M×N), thus obtaining physical property models;normalizing the physical property models, thus obtaining normalizedphysical property models ND_(data) (L,M,N) and NR_(data) (L,M,N) whichare then represented in a scatter plot of physical properties;converting distribution positions of the ND_(data) (L,M,N) and NR_(data)(L,M,N) on the scatter plot of physical properties into type angle (TA)using a four-quadrant inverse tangent and into type intensity (TI) usinga distance from an origin; determining a minimum of TI values, whichclassifies two or more classes on a scatter plot of the TA and TIvalues, to be a threshold; and analyzing the 3D geological structurebased on a local extreme point and/or points in the scatter plot for theTA and TI using the threshold, wherein analyzing the 3D geologicalstructure is performed by determining a rock type, and a structuredepending on physical properties and geological positions of theclasses.
 2. The method of claim 1, wherein the normalizing is determinedby Equation 1 below:ND=[D−Mean(D)]/[(Max(D)−Min(D))/2]  <Equation 1> wherein D is a physicalproperty model D_(data)(L,M,N) or R_(data)(L,M,N), and ND is anormalized physical property model.
 3. The method of claim 2, furthercomprising comparison which is performed in such a manner that whetheror not the two or more classes overlap with each other is determinedbased on the determined threshold, between determining the threshold andanalyzing the 3D geological structure.
 4. The method of claim 3, whereinwhen the two or more classes overlap with each other under conditions ofbeing equal to or more than the threshold, determining the threshold isrepeated so that the threshold is re-determined.
 5. The method of claim1, further comprising, before estimating the density value on the commoncoordinates: surveying gravity at a position at which the 3D geologicalstructure is to be analyzed, thus obtaining gravity data; correcting thegravity data, thus extracting a Bouguer anomaly; and performing 3Ddensity inversion using the Bouguer anomaly, thus obtaining the densitymodel.
 6. The method of claim 5, further comprising, before estimatingthe resistivity value on the common coordinates: surveyingmagnetotelluric (MT) at a position at which the 3D geological structureis to be analyzed, thus obtaining MT data; correcting the MT data basedon a remote reference, thus acquiring extracted MT data; and performing3D MT inversion using the extracted MT data, thus obtaining theresistivity model.
 7. The method of claim 6, wherein the density modeland the resistivity model on the common coordinates are estimated usingresults analyzed via inversion after the common coordinates have beenset in the density model and the resistivity model and using ordinarykriging upon estimating the density value and the resistivity value onthe common coordinates.
 8. A method of analyzing a three-dimensional(3D) geological structure, comprising: estimating a first model valueand a second model value on common coordinates L, M, and N to calculatea first model D_(data) (L,M,N) and a second model R_(data) (L,M,N) on asame 3D grid (L×M×N) having different physical properties, thusobtaining physical property models; normalizing the physical propertymodels, thus obtaining normalized physical property models ND_(data)(L,M,N) and NR_(data) (L,M,N) which are then represented in a scatterplot of physical properties; converting distribution positions of theND_(data) (L,M,N) and NR_(data) (L,M,N) on the scatter plot of physicalproperties into TA using a four-quadrant inverse tangent and into TIusing a distance from an origin; determining a minimum of TI values,which classifies two or more classes on a scatter plot of the TA and TIvalues, to be a threshold; and analyzing the 3D geological structurebased on a local extreme point and/or points in the scatter plot for theTA and TI using the threshold, wherein analyzing the 3D geologicalstructure is performed by determining a rock type, and a structuredepending on physical properties and geological positions of theclasses.
 9. The method of claim 8, wherein the normalizing is determinedby Equation 1 below:ND=[D−Mean(D)]/[(Max(D)−Min(D))/2]  <Equation 1> wherein D is a physicalproperty model D_(data)(L,M,N) or R_(data)(L,M,N), and ND is anormalized physical property model.
 10. A method of analyzing athree-dimensional (3D) geological structure, comprising: estimatingmodel values on common coordinates L, M, and N to calculate two or moremodels on a same 3D grid (L×M×N) having different physical properties,thus obtaining physical property models; normalizing the physicalproperty models, thus obtaining two or more normalized physical propertymodels which are then represented in a scatter plot of physicalproperties; converting distribution positions of the normalized physicalproperty models on the scatter plot of physical properties into TA usinga four-quadrant inverse tangent and into TI using a distance from anorigin; determining a minimum of TI values, which classifies two or moreclasses on a scatter plot of the TA and TI values, to be a threshold;and analyzing the 3D geological structure based on a local extreme pointand/or points in the scatter plot for the TA and TI using the threshold,wherein analyzing the 3D geological structure is performed bydetermining a rock type, and a structure depending on physicalproperties and geological positions of the classes.
 11. The method ofclaim 10, wherein the normalizing is determined by Equation 1 below:ND=[D−Mean(D)]/[(Max(D)−Min(D))/2]  <Equation 1> wherein D is a physicalproperty model, and ND is a normalized physical property model.